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ABSTRACT 

Cold collapse of a cluster composed of small identical clumps, each of which is in virial 
equilibrium, is considered. Since the clumps have no relative motion with respect to 
each other initially, the cluster collapses by its gravity. At the first collapse of the 
cluster, most of the clumps are destroyed, but some survive. In order to find the 
condition for the clumps to survive, we made systematic study in two-parameter space: 
the number of the clumps Nc and the size of the clump r^ . We obtained the condition, 
Nc S> 1 and Uk > 1, where rik is related to r^ and the initial radius of the cluster 

Rini through the relation Rmi/rv = 2Nc'' '' . A simple analytic argument supports 
the numerical result. This n^ corresponds to the index of the power spectrum of the 
density fluctuation in the cosmological hierarchical clustering, and thus our result 
may suggest that in the systems smaller than 2/(ri/i^)Mpc, the first violent collapse 
is strong enough to sweep away all substructures which exist before the collapse. 
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1 INTRODUCTION 

Hierarchical structures are seen in various astrophysical ob- 
jects in the universe, such as globular clusters in a galaxy, 
or galaxies in a cluster of galaxies. Such a subsystem in a 
virialised system forms a distinct object which is in virial 
equilibrium by itself. How such subsystems form and evolve 
in a larger scale environment is a fundamental problem of 
astronomy. 

This is because, in general, mutual gravitational inter- 
action tends to destroy the subsystems. Once a composite 
system of a large halo and many subsystems forms, there are 
several destructive processes, such as tidal heating caused by 



For clusters of galaxies, the main effective processes to 
disrupt galaxies are encounters between the galaxies and 
the tidal heating by the enco unter with the cluster core. 
Richstone & Malumath (1983) calculated the evolution of 
clusters of galaxies, taking into account the tidal stripping 
and mergers of the galaxies with a Monte-Carlo method and 
showed that most galaxies could survive until a Hubble time, 
despite large mass loss. By using the numerically determined 
cros s secti on for a merger of two identical galaxies, Makino & 
Hut( 1997| ) estimated the merger rate of galaxies in a cluster, 
and concluded that a galaxy experiences mergers, on aver- 
age, once every 100 Gyr. On the other hand, some numeri- 



and a core of Llie halu( gpilzcr 1908 ). For Llie easy o f }j,lobular 



cal s imulations (Barnes 1989; Funato, Makino & Ebisuzaki 



^los^ ei]counters amon g suDsystcm s oj Between a suDsystcm ^gggp ^^^^ that mergers are very effective in a cluster with 



clusters in a galaxy, bulge shocking, disc shocking ( Ostrikerj 
B pitzer fc Chevalier 1974 ) ^nd stellar evaporation( ppitzer & 
Hart 1971) are also effective. 

The rate of the disruption of subsystems depends on 
the properties of the subsystems and the environment. In 
the Galaxy, for example, a cluster should have a size and 
mass in a ce rtain narrow range in or der to survive for a 
Hubble time(|Fall fc Rccs 1977|; |Caputi 



Castellani 1984 



Chernoff fc Shapiro 1987] ; pkazaki fc Tosa 1995| ). Moreover 
the surviv al factor depend s on the distance from the Galac- 
tic centre( Vesperini 1997 ). More than a half of the present 
globular clust ers in the Galaxy will be destroyed in the next 



Hubble tim e (Aguilar, Hut 
triker 19971). 



Ostriker 198f 



Gnedin & Os- 



a relatively small number of galaxies {^ 100); consequently 
the galaxies will disappear rapidly in a Hubble time. 

It might be noted that these studies considered the evo- 
lution of the virialised subsystems in a completely virialised 
cluster. This is not the case, however, in a hierarchical clus- 
tering scenario, since smaller systems form earlier so that 
galaxies form earlier than clusters of galaxies. Hence the 
galaxies in a cluster experienced the early violent phase of 
the cold collapse of the cluster. Also for globular clusters, 
there is a theory that gl obular clusters fo rms during the 
collapse of a protogalaxy (Fall & Rees 1985). In both cases, 
substantial subsystems should have destroyed in the dynam- 
ically evolving phase of the parent system. 

Though destructive processes in subsystems should be 
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more effective in a collapsing system than in a virialised one, 
only a few studies have been made for the former and we 
have not yet arrived at a fundamental understanding. This is 
partly because analytic approaches are inapplicable to dy- 
namically evolving systems. In addition, this is a difficult 
problem also for numerical calculations, since the existence 
of compact subsystems requires high resolution in mass and 
time, which means that we need a large number of particles 
and a small time step to cope with such systems. As a re- 
sult it takes huge computational time. A pionee ring study 
of this problem was made by van Albada ( 1982 ). He stud- 
ied cold collapses numerically with a wide variety of initial 
conditions in order to make clear the origin of the observed 
light distribution of elliptical galaxies. In his simulations, 
however, clumps are so large that their sizes are compara- 
ble to the mean separation of the clumps. Hence all clumps 
would be disrupted by the cold collapse. Later, the impor- 
tance of the present issue was recognised in cosmological 
simulations of structure formation. In the numerical simu- 
l ations, substructures are swept away in clusters of galaxies 



Frenk et al. 198J)- It is called the "over- 



( [Whiteet al. 1987 ; 

merging problem". Moore, Katz & Lake (199(3) showed that 
current numerical resolution is sufficient to prevent dissolu- 
tion by two-body relaxation, and that physical effects, such 
as tidal heating, are important. They also argued that the 
size and the density of the innermost core of the subsystems 
are crucial factors for their survival, and the current numeri- 
cal resolution may not be enough to resolve the overmerging 
problem. 

In this paper we investigate the fundamental destruc- 
tion mechanisms of the subsystems, especially in the circum- 
stance that the environment is dynamically evolving, using 
both numerical and analytical approaches. We employ a sim- 
ple model, in which the system initially consists of identical 
virialised clumps. The clumps have no relative motion, so 
that the system collapses by gravity. This model is charac- 
terised by only two parameters: the number of clumps in 
a cluster and the size of the clumps relative to that of the 
cluster. We make a systematic numerical survey in the two- 
parameter space. 

In Section 2, we explain the initial conditions of our 
simulations. The results of the simulations are presented in 
Section 3. A simple analytic estimate of the condition for 
the clumps to survive is given in Section 4. Final section is 
devoted to the summary. 



2 INITIAL CONDITIONS AND MODEL 
PARAMETERS 

Our main aim is to find the condition for the survival of sub- 
structures during a cold collapse. Qualitatively, we expect 
that smaller clumps are more likely to survive because their 
binding energies are large. In addition, since the collisional 
effect between stars in a stellar system becomes weakened 
as the number of the stars increases, we also expect more 
clumps to survive as the number of clumps increases. 

The dependence of the fraction of survival on the num- 
ber of the initial clumps is rather uncertain. In order to find 
a quantitative condition, therefore we need to survey all the 
parameter space systematically by adopting simple initial 
conditions. 



The initial models are characterised only by two param- 
eters; one is the size of clumps and the other is their number. 
All clumps are identical in mass and size, and all particles 
belong to one of the clumps. 

Each clump is realised as an equilibrium Plummer 



model (Binney & Tremaine 1987); the density distribution 



p(r) is given by 
3mi 



P = 



47rrf 



l + (-) 



-5/2 



(1) 



where r is the radius, ri is the core radius, and mi is the 
mass of the clump. The potential and distribution function 
of a particle at the position, x, with the velocity, v, are 
given, respectively, by 

Gmi 



<I>(r-) = - 



and 



^Jrl+■. 



(2) 



f{^,^)=K\e\"\ (3) 

where e is the specific energy of a particle, e = ^«'^-|-$(r), G 
is the gravitational constant, and is' is a numerical constant. 
Since the smallest structure and shortest time scale are 
determined by those of the clump, we take the mass, size, 
and dynamical time of the clump as the units of computa- 
tion. Specifically, G = m\ = r„ = 1, where r„ is the virial 
radius defined by energy of the clump. 



E^- 



Gml 

4r„ 



(4) 



In these units, the half mass radius is ru = 0.769, and the 
crossing time at rn is 0.834. 

A cluster is composed of Nc clumps, which are dis- 
tributed uniformly within a sphere. The centre of mass of 
each clump is located at a point which is randomly selected 
from the sphere and has no bulk velocity initially. Since the 
size of the clumps is fixed, the initial radius of the cluster 
Rini is scaled by 



Ri; 



= 2N, 



(7ifc+5)/6 



(5) 



where n^ is a parameter. This relation corresponds to the 
ratio between the radii of virialised objects in the Einstein- 
de Sitter universe with the initia l pow er spectrum such as 
P{k) (X k"'' [e.g., §26 of Peebles ( |l980| )]. The larger rik be- 
comes the smaller relative size of the clumps. In particular, 
when Uk = —3 the mean density of the cluster is about the 
same as that of a clump. 

In our simulation, the number of particles in a clump is 
fixed to be 512. With this number, the clump is stationary 
at least for T ^ 50. In the longest run two-body relaxation 
may be effective, b ut str uctural change is very small. Fur- 
ther, Moore et al. (1996) assure that for A'^ ^ 100 numeri- 
cal effects on disruption of clumps are negligible compared 
with physical effects such as tidal heating. For confirma- 
tion, we made test runs with a half (A^ = 256) and a double 
{N = 1024) particles and found the difference in the number 
of surviving clumps is within a factor of 15%. 

We examine the cases with four different values of Uk 
: Uk — —3, —2, —1, 0, and for each Uk, we assign eight 
different values of N^ : Nc = 2, 4, 8, 16, 32, 64, 128, 256, 
except for Uk — 0. Thus we have 31 different combinations 
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of Tife and Nc- Further we prepare two different realisations 
of random distribution of clumps, which are referred to as 
RUNl and RUN2. In Table 1, we summarised the various 
parameters adopted in the present simulations. The first and 
second columns are n^ and Nc, respectively. The third, the 
fourth, and the fifth columns are filling factor of clumps in 
the whole system, the initial radius and virial ratio of the 
system, respectively. These three parameters are determined 
according to Uk and N^. 

In our numerical code, gravity is calculated by direct 
sum mation of contributio n of all particles using GRAPE- 
3A ( Okumura et al. 1993 ). The softening length, e, and the 
time step of integration, St, are set to be e = 1/32 and 
St = 1/64, respectively. The total energy of the system is 
conserved within an accuracy of a few hundredth of 1% on 
average and at maximum 1%. 



3 NUMERICAL STUDY FOR SURVIVAL 
CONDITION 

The cluster starts collapsing in a similar way to the case of 
a cold collapse without substructure. Since the distribution 
of clumps in the cluster is spherical and uniform, all the 
clumps collapse simultaneously at T = Tcoii- In this paper, 
we refer to this epoch as the first collapse. The values of Tcoii 
for different parameters are listed in Table 1. Figs, nl and Q 
show two examples of evolution : RUNl with nt ~ —3 and 
Nc = 64 (example 1), and RUNl with n^ = and N^ = 64 
(example 2) . We follow the evolution until T = 2 Tcoii , be- 
cause we are concerned only with the first violent stage. 
After T = 2Tcoii, the system no longer experiences violent 
change in the potential, thus the processes of destruction of 
clumps are the same a s those which have alreadv be e n stud- 
ied b y many authors ( Richstone fc Malumath 1983 ;_ Barnei 



198£; [Funato, Makino fc Ebisuzaki 1993t [Makino fc Hul 



1997). The final state of our simulation could be the initial 
state of the above studies. At the end of our simulations, we 
count the number of the clumps which survive through the 
first violent stage. 

In the example 1, all clumps are destroyed at the first 
collapse and merge into a large smooth halo. In the example 
2, in contrast, many clumps still have their identities even 
at T = 2 Tcou ■ (Note that most of the survivors are in the 
central core, but they cannot be resolved in Fig. bl) 

Kinematical diagrams are much more useful to identify 
the surviving clumps rather than density distribution. If a 
system is in a dynamical equilibrium, the motion of a par- 
ticle in the phase space is periodic so that the trajectory 
should look like an oval or a spindle in a position-velocity 
plane. Fig. shows the evolution of the example 1, in the 
(y-Vy) plane, as a typical example. We examined all three 
planes of the (x-v^c), (y-Vy), and (z-Vz) plane and pick up 
the clumps with the spindle structure in the all three plane. 

At the beginning (see Fig. wk) each clump forms an 
individual spindle, whereas in the final state (see Fig. Bp) 
many of the spindle structures are destroyed during the first 
collapse, and a large merged halo is formed at the centre. 
There exist in the latter two prominent features: the central 
big spindle structure, which is due to the merged halo, and 
the rectilinear wings by escaping particles from the merged 
halo. Besides, there are small spindle structures superim- 



posed on the main structures. These are due to the surviv- 
ing clumps. Fig. H demonstrates the different properties of 
the surviving and disrupted clumps. Two clumps are cho- 
sen as representatives. Clump #16 is destroyed in the first 
collapse and spreads over the merged halo afterwards. Con- 
versely, Clump #17 still retains the equilibrium structure 
a the end of calculation. We count the number of all the 
surviving clumps after the calculation, and list them in the 
column labelled A'survivc of Table 1. 

There are some clear tendencies in the results: first, ir- 
respective of the index Uk, the cluster with a small number 
of clumps {Nc < 16) yields no surviving clumps. Second, 
for Nc > 32 some clumps may survive. The number of the 
survivors depends on nt- For rife = —3, for example, all the 
clumps are destroyed except for the case with Nc = 128. For 
the larger Uk and the larger A''c,the more clumps survive. 

One thing which we should note when considering the 
mechanism for the disruption of the clumps is mass loss or 
escapers. It is known from the studies of the cold collapses 
with uniform initial density distribution without clumps, 
that models with a small initial virial ratio are liable to 
suffer substantial mass loss. In the cases with clumpy initial 
distribution, we observed even more mass loss. Among the 
escaping particles, some clumps escape to the infinity and 
some not. It is thus better to discriminate unbound sur- 
vivors from bound ones in the central big halo, since the 
subsequent evolution is distinct for those two. As unbound 
survivors, we pick up those which lie in the rectilinear wings 
in the kinematical diagram. In Fig. H, for example, four such 
clumps can be seen in the wing extending to the negative 
y direction. They are escaping from the central halo with a 
constant velocity. Some other clumps are escaping in the x 
and z directions, but still many survivors are found to orbit 
in the central halo. The numbers of the bound survivors are 
listed in the column labelled A'bound of Table 1. The num- 
ber of the bound survivors increases as rit and Nc increase. 
However, the largest difference between the numbers of the 
bound and all survivors can be seen in Uk = —2. Only few 
survivors are found to be bound, while there are many es- 
caping survivors. For rife = —1 and 0, the number of bound 
survivors increases as A'^c increases as long as Nc ^ 64. 

Thus we now derive the conditions on the initial pa- 
rameters for the cluster to contain distinct clumps when it 
becomes virialised: 



Uk > -1 and Nc Z 64. 



4 SIMPLE ANALYTIC MODEL 



(6) 



In order to understand the physical mechanism underlying 
disruption of clumps, we construct a simple analytic model. 
Suppose that a cold cluster with Nc clumps start to col- 
lapse by its own gravity, same situation as those described 
in the previous section. We surmise that the most effective 
mechanism for destroying clumps is the tidal heating by the 
gravity of the whole cluster, especially at the maximal col- 
lapse, where the radius of the system becomes minimum. Let 
m and r„ be the mass and the the virial radius of a clump, 
respectively, M — Nciri be the total mass of the cluster, and 
Rmin be the minimum radius of the cluster at the maximal 
collapse. The clumps move in the gravitational potential of 
the cluster with a typical distance Rmin and with a typical 
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velocity V ~ {GM/R-minY ■ Then the change in the en- 
ergy of a clump is given by the ti dal h eating formula [e.g., 
eq. (7-55) of Binney & Tremaine (|l98^)], 



A£i 



G^Arm 2 _ GM 

mtn mm 



mr,, = 4 






where E\ is the energy of the clump given by 
Grr? 



£i = — 



4r„ 



\Ei\, (7) 



(8) 



In the above formula, the minimum radius Rmin is an 
unknown paramet er. He re we employ the result of Aarseth, 
Lin & Papaloizou ( 1988 ) ; in the cold collapse of a sphere with 



a uniform distribution of A^ particles, the system contracts 
to the minimum radius Rmin, depending on the number of 
particles. 



^min — J' ^ini, 



(9) 



where Rini is the initial radius of the sphere. Since this min- 
imum radius results from the statistical A^" fluctuation, 
we can apply the result to our system by substituting the 
number of clumps Ac for A in equation (bl). Moreover, using 
equation (pi) as a relation between Rini and r^, the relative 
change in the energy of the clump can be expressed only in 
terms of Ac and Wfc: 



ASi 



= 4Ac A, 



-1/3 



2a; 



("fc+5)/6] 



2 " 



("fc+l)/2 



(10) 



If rifc -I- 1 is positive, the relative change in the energy de- 
creases as Ac increases. When Ai?i/|_Ei| is sufficiently small, 
a large fraction of clumps are expected to survive the first 
collapse. 

Surprisingly, this simple estimate is in excellent agree- 
ment with the results of the numerical simulations, espe- 
cially for the number of the bound survivors. Therefore, we 
conclude that the disruption of the clumps is mainly due 
to the tidal shock at the first collapse, and that the condi- 
tions on the initial parameter for the cluster to sustain the 
substructure until reaching the virial equilibrium are 



rik 



> -1 and Ac > 1. 



(11) 



Note the estimate (jlOl) is valid for those clumps which 
collide with each other simultaneously at the first col- 
lapse. Owing to the statistical fluctuations, there exist some 
clumps which collapse belatedly. Those clumps are acceler- 
ated to a larger velocity, and some may gain enough energy 
to escape to the infinity. When the delayed clumps reach 
their pericentre, however, the cluster as a whole has already 
proceeded into the reexpansion phase, and thus the gravi- 
tational force is weaker than that at the maximal collapse. 
Therefore, the escaping clumps are more likely to survive. 
In fact, we found in the simulations a tendency that the 
escaping survivors are initially located in the outer part of 
the cluster. The bound survivors have a similar tendency as 
well. 



index Uk of the spectrum of initial density fluctuation as 
two parameters characterising the initial distribution, the 
condition is expressed as 



Tifc > -1 and Ac > 1. 



(12) 



As an application of this result, let us consider, for exam- 
ple, the CDM model of cosmological clustering. It is known 
that in the CDM model the power spectrum of the density 
fluctuation [eq. (25.22) of Peebles(1993)] is given by 

P{k) oc ^, (13) 

[1 + 8fc/t7/i2 + 4.7fcV(!r2/i2)2]2 ' 

and the gradient of the power spectrum is larger than 
Uk — —1 for the fluctuation with the scale larger than 
2/(il/i^)Mpc. Therefore our result may suggest that in a 
system smaller than 2/{Q,h'^)Mpc, the flrst violent collapse 
is strong enough to sweep away all substructures which exist 
before the collapse. 

Since our simple initial mass model has, in fact, no spec- 
trum in the distribution of masses of subsystems, it is more 
useful to express this condition in terms of other parame- 
ters. Using radii of the cluster and the clumps, for example, 
we derive 

^ini 



Ac > 1 and 



^2Ac/^ 



(14) 



For a system with many clumps, the flUing factor C, 
which is the inverse of the fraction of the volume occupied 
by the clumps, is widely used to describe the dumpiness. In 
our case, C = Ri„i/{Ncry), so that the condition is 



Ac > 1 and G^ 8Ac 



(15) 



Since our numerical model is quite idealised, we need to 
incorporate many factors for more realistic discussion. First, 
as the model of each clump, we examine only the Plummer 
sphere given by equation (hi). This model has a core whose 
radius is about a half of the half mass radius, r^. On the 
other hand, the observed globular clusters or elliptical galax- 
ies ha ve core radii of about 1/10 or 1/100 of the half mass 
radii (Binney & Tremaine 1987). Such smaller and denser 
cores are expected to survive even when haloes or envelopes 
are stripped. Second factor is rotation or asymmetry. Both 
make the flrst collapse of the cluster less dense, and thus 
more clumps should survive after the collapse of the cluster. 
All these factors weaken the destructive effect. Therefore, 
the condition, presented in this paper, for substructure to 
survive should be sufflcient condition. 
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5 CONCLUSIONS AND DISCUSSION 

We have obtained the condition for substructures to sur- 
vive after the whole system has collapsed to reach a virial 
equilibrium. With the number of clumps, Ac, and the power 
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Table 1. List of parameters and results of our numerical simulations 



nk 


No 


C 


Rini 


virial ratio 


TcoU 


N 

-'^ survive 


-'^bound 


-3 


2 


8. 


2.52 


0.73 


3.6 


0/0 


0/0 


-3 


4 


8. 


3.17 


0.57 


4.8 


0/0 


0/0 


-3 


8 


8. 


4.00 


0.25 


4.0 


0/0 


0/0 


-3 


16 


8. 


5.04 


0.24 


4.0 


0/0 


0/0 


-3 


32 


8. 


6.35 


0.14 


3.6 


0/0 


0/0 


-3 


64 


8. 


8.00 


0.099 


3.6 


0/0 


0/0 


-3 


128 


8. 


10.1 


0.063 


3.4 


2/2 


0/0 


-3 


256 


8. 


12.7 


0.040 


3.3 


0/6 


0/0 


-2 


2 


11.3 


2.83 


0.75 


4.2 


0/0 


0/0 


-2 


4 


16. 


4.00 


0.62 


6.6 


0/0 


0/0 


-2 


8 


22.6 


5.66 


0.43 


6.5 


0/0 


0/0 


-2 


16 


32. 


8.00 


0.32 


7.9 


2/0 


0/0 


-2 


32 


45.3 


11.3 


0.23 


8.4 


0/2 


0/0 


-2 


64 


64. 


16.0 


0.18 


10.0 


4/2 


1/0 


-2 


128 


90.5 


22.6 


0.13 


11.2 


12/10 


2/0 


-2 


256 


128 


32.0 


0.095 


13.2 


25/30 


0/4 


-1 


2 


16. 


3.17 


0.77 


4.8 


0/0 


0/0 


-1 


4 


32. 


5.04 


0.67 


9.2 


0/0 


0/0 


-1 


8 


64. 


8.00 


0.51 


10.9 


0/0 


0/0 


-1 


16 


128. 


12.7 


0.43 


15.0 


0/1 


0/0 


-1 


32 


256. 


20.2 


0.35 


19.5 


3/2 


0/0 


-1 


64 


512. 


32.0 


0.30 


28.0 


18/19 


8/10 


-1 


128 


1024. 


50.8 


0.25 


37.8 


37/35 


13/17 


-1 


256 


2048. 


80.6 


0.21 


52.0 


80/67 


23/24 





2 


22.627 


3.56 


0.79 


5.6 


0/0 


0/0 





4 


64 


6.35 


0.72 


12.9 


0/0 


0/0 





8 


181.01 


11.3 


0.60 


18.1 


0/0 


0/0 





16 


512 


20.2 


0.55 


30.5 


1/0 


1/0 





32 


1448.1 


35.9 


0.49 


46.0 


5/4 


4/2 





64 


4096 


64.0 


0.46 


79.2 


24/23 


18/11 





128 


11585. 


114.0 


0.43 


126.7 


70/54 


48/35 



Each column means: (l)nfe the index of the power spectrum, {2)Nc the number of clumps, (3) C filling factor, (4) Rini 
initial radius of the cluster (5) the initial virial ratio, (6) T^-aii the time of the first collapse, (7)Afsurvive the number of 
clumps which have a bound core at the final state for RUNl and RUN2, (8) Af^ound the number of surviving clumps which 
is gravitationally bound in the merged halo for RUNl and RUN2. 
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Figure 1. Time evolution of the projected distribution in the (j/, z) plane of the system with nj; 
time is indicated at the upper-right corner of each figure. 
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Figure 2. Same as Fig. hi but the case with nj. = and Nc = 64 (example 2). 
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Figure 3. Kinematical diagram of the systems with nj. = and Nc = 64 in the (y-Vy) plane (a)at the beginning (T = 0), and (b) at 
the end of the calculation (T = 160). 
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Figure 4. Kinematical diagram of the same systems as in Fig. pi but showing only two clumps as an example. (a)At the beginning 
(T = 0), and (b) at the end of calculation (T = 160). The structure of Clump #16 does not change much while Clump #17 is destroyed 
and spreads over the merged halo. 
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